Draft version August 21, 2012 

Preprint typeset using I^-T^X style emulatcapj v. 5/2/11 



(N 

o 

(N 
< 

o 



Of 

6 

o3 



> 

(N 
O 

00 

o 

(N 



13 



SCATTERING OF THE F-MODE BY SMALL MAGNETIC FLUX ELEMENTS FROM OBSERVATIONS AND 

NUMERICAL SIMULATIONS 

T. Felipe 1 , D. Braun 1 , A. Crouch 1 , A. Birch 1,2 

Draft version August 21, 2012 

ABSTRACT 

The scattering of /—modes by magnetic tubes is analyzed using three-dimensional numerical sim- 
ulations. An /—mode wave packet is propagated through a solar atmosphere embedded with three 
different flux tube models which differ in radius and total magnetic flux. A quiet Sun simulation 
without a tube present is also performed as a reference. Waves are excited inside the flux tube and 
propagate along the field lines, and jacket modes are generated in the surroundings of the flux tube, 
carrying 40% as much energy as the tube modes. The resulting scattered wave is mainly an /—mode 
composed of a mixture of m = and m = ±1 modes. The amplitude of the scattered wave approx- 
imately scales with the magnetic flux. A small amount of power is scattered into the p\ — mode. We 
have evaluated the absorption and phase shift from a Fourier-Hankel decomposition of the photo- 
spheric vertical velocities. They are compared with the results obtained from the emsemble average of 
3400 small magnetic elements observed in high-resolution MDI Doppler datacubes. The comparison 
shows that the observed dependence of the phase shift with wavenumber can be matched reasonably 
well with the simulated flux tube model. The observed variation of the phase-shifts with the azimuthal 
order m appears to depend on details of the ensemble averaging, including possible motions of the 
magnetic elements and asymmetrically shaped elements. 
Subject headings: MHD; Sun: oscillations 



1. INTRODUCTION 

The work by iBraun et al.1 ()1988l ) has shown that 
sunspots can absorb up to half the power of incident 
p— modes. At the same time, part of the wave flux is 
scattered. The p— mode absorption and scattering phase 
shifts depend on the frequency, degree, radial order 
and azimuthal ord er of the incident mode (Bogda n et al.l 
119931 lBraunlll995l ). and also on the magnetic structure, 
making the study of the scattering a promising way to 
infer the subsurface structure of sunspots and other mag- 
netic features. In these studies the authors used Hankcl 
analysis, a method which decomposes the p— modes into 
inward and outward propagating waves in annuli sur- 
rounding the sunspot. 

Several mechanisms have been propos ed to explain the 
obser ved absorption, mode conversion (ICallv fc Bogdanl 
119931 ) being the principal candidate. When an acous- 
tic wave encounters a magnetic field concentration, it is 
split into fast and slow modes. This mode transforma- 
tion occurs at the height where the Alfven and sound 
velocities are comparable, since at that layer the dis- 
tinction between the modes is small. Below this height 
the sound speed is higher than the Alfven speed and 
the modes are effectively decoupled. The fast mode is 
an acoustic-like wave, while the slow mode is similar to 
an Alfven wave and propagates along field lines remov- 
ing energy from the acoustic wave. The observed ab- 
sorption of /— modes can be acco unted for by a vertical 
magnetic field (|Callv et al.l 119941) . but the p— mode ab- 
sorption obtained from this model is insufficient. How- 
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ever, the presence of inclined magnetic field produces sig- 
nific ant increases in absorp tion with a peak at around 
30° (jCrouch fe Callvl 120031) which are consistent with 
observed values ICallv et al.l[2003l ). Other mechanisms 
may also play a role. One of them is resonant absorp- 
tion (|Hollweelll988l : IRosentha]lll992T) . It may occur when 
the flux tube has a smooth variation of the magnetic 
field rather than discontinuous, and represents absorp- 
tion of wave energy by the transition layer when the inci- 
dent acoustic waves resonantly excite MHD waves in the 
magnetic structure. However, the amount of absorption 
achieved by this mechan ism cannot exp lain the observed 
loss of p— mode power. iD'Silval (|1994[ ) points out that 
apart from this absorption produced by the dissipation 
of the p— modes in resonant layers and mode conversion, 
mode mixing also takes places. In this process the power 
of an incident p— mode mode with a certain frequency 
and degree I can be dispersed into an outgoing p— mode 
with the same frequency, but different degree. 

Flux tubes are a key feature to understand solar mag- 
netic activity. They are spread all over the solar surface 
and couple different layers of the solar atmosphere. It has 
been proposed that magnetic flux tubes can act as wave 
guides, being one of the possible sources which supply en- 
ergy to the upper lay ers to account for the chromospheric 
and coronal heating (|Jefferies et al.l 12006) . The interac- 
tion of p— modes with thin flux tubes excites tube waves, 
including sausage waves, which are axisymmetric, longi- 
tudinal waves driven by variations in the total pressure, 
and kink waves, whose restoring force is magnetic tension 
and buoyancy, producing transversal oscillations. These 
waves propagate upward or downward and extract energy 
from the p— modes of the acoustic cavity dBogdan et all 
119961 : iHindman fe Jainll2008l : I Jain et alJl2009f ). The kink 
mode is driven by the distortion of the tube produced by 
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the harmonic flow field of the p— modes, while the sausage 
wave is excited by th e pressure pertubati ons associated 
to the acoustic waves (|Bogdan et aljfl996t) . These mech- 
anisms arc different from mode conversion, discussed in 
the previous paragraph, where fast and slow magnetoa- 
coustic waves exchange energy due to their strong cou- 
pling in the region wher e the sound and Alfven speeds are 
comparable. Recently, Daiffalla h et al.l (|2011l) used nu- 
merical simulations to study the scattering of an /—mode 
by vertical flux tubes of different sizes, finding that the 
scattering by tubes with small radius is dominated by 
the kink mode, while the sausage mode is dominant for 
lar ge tubes. Thi s resu lt coincides with the earlier work 
bv lBogdan et al.l (|1996l ). who studied analytically the na- 
ture of wave interactions with thin flux tubes (|Spruitl 
1 1 9 8 lh and found that the kink mode is the dominant tube 
wave. The thin flux tube approximation assumes that 
the diameter of the tube is smaller in comparison to the 
pressure scale height, and t hus the horizontal vari ations 
inside it can be neglected. lHanasoge et al.l (|2008l ) eval- 
uated the scattering matrix associated with an /—mode 
that interacts with a thin flux tube in a stratified atmo- 
sphere, focusing on the kink mode excited in the mag- 
netic tube. They found that most of the scattered wave 
corresponds to an /—mode with amplitude of 1 .17% and 
with a phase shift of abound 50° relative to the incident 
wave, overstimatin g the observed value by a factor of 8.8 
(IDuvall et al.ll2006f). On th e other hand, the recent work 
bv iHindman fc Jainl ()2012h analyzed the axisymmetric 
scattering of p— modes, mediated through the excitation 
of sausage waves on the flux tube, instead of the kink 
mode. They obtained a small absorption due to the poor 
coupling between the /—mode and the sausage mode for 
thin flux tubes, as pointed out by some of the works pre- 
viously described in this paragraph. 

Although these theoretical works have provided the 
first predictions about the modification of the solar wave 
field produced by flux tubes, as far as our concern no 
attempt has been made to obscrvationally measure the 
detailed properties of the scattering produced by these 
small magnetic elements, with the exception of the es- 
timates of a mplitude and phase of monopole and dipolc 
scattering bv IDuvall et a l. (2006). One of the objectives 
of this work is to present the measurement of the phase 
shift and its variations with the azimuthal orders m and 
degree L. This data is a fundamental input to confront 
with the theory. On the other hand, most theoretical 
studies of this topic have been based on an analytical 
development. All of them have been restricted by some 
limitations, including the use of the thin flux tube ap- 
proximation, the lack of the gravitational stratification, 
the analysis of a polytropc instead of a realistic solar 
atmosphere, or some constrains in the process that me- 
diates the scattering. Numerical simulations are a more 
versatile approach and allow us to study more general 
situations. 

In this work we study the scattering of an /—mode 
by flux tubes of different radius and magnetic flux using 
numerical simulations. As discussed previously, it is well 
known that thin flux tubes support sausage and kink 
modes (|Roberts fc Webbl [l978t ISpruitl [l98l . On the 
other hand, in unstratified atmospheres permeated by 
homogeneous magnetic fields one would expect the prop- 
agation of pure fast and slow magnetoacoustic waves and 



Alfven mode. In atmosphere stratified by gravity (for ex- 
ample) the fast, slow, and Alfven waves are coupled in 
general and this distinction between modes no longer ap- 
plies, but even in these cases it is useful to refer to this 
simple picture to discuss the properties of fast and slow 
magnctoacoustic-like waves (in regions where the sound 
and Alfven speed differ greatly). In the case of thick 
flux tube models, representative of a sunspot, for exam- 
ple, acoustic waves can be converted into these modes 
by means of mode conversion. The flux tube models pre- 
sented in this paper correspond to an intermediate case 
between these two extremes. We expect a smooth tran- 
sition (with increasing radius) from mostly excitation of 
the kink and sausage modes at small radius, to excitation 
of waves that look like the fast and slow magnetoacoustic 
waves of a thick flux tube (at large radius) . However, in 
this study we made no attempt to distinguish between 
neither the different wave modes which are present nor 
the mechanisms that generates them. Instead, we will 
refer to "tube modes" or "tube waves excitation" indis- 
tinctly. We aim to carry out a direct comparison between 
the numerical and observational results by performing a 
Hankel analysis of the data obtained from the interaction 
of an /—mode with flux tubes. The organization of the 
paper is as follows. In Section [2] we introduce the obser- 
vations used in this study. Section [3] briefly describes the 
numerical code and the set up of the simulations. In Sec- 
tions |4] and [5j we present the tube mode excitation, the 
scattering, and the jacket modes produced in these sim- 
ulations. The results of the Hankel analysis are shown in 
Section^ including a comparison with observations, and 
finally we conclude with a summary of our calculations 
and a discussion of their applicability to understanding 
observations. 

2. OBSERVATIONS 

The motivation of this paper is to understand how 
small magnetic elements affect /-mode wave propaga- 
tion. The comparison of observations of magnetic el- 
ements on the Sun with equivalent measurements ob- 
tained from numerical simulations can help infer the 
properties of the solar magnetic elements. These com- 
parisons can also be used to assess the limitations of, 
and potentially improve, the observational techniques. 

Observations of solar magnetic elements, using Dopp- 
lcrgra ms obtained from th e Michelson Doppler Imager 
(MDI; iScherrer et al.l Il995h onboard the Solar and He- 
liospheric Observatory (SOHO), were analyzed using the 
Fourier-Hankel spe ctral decomp osition method, as de- 
scribed in detail by Braun (1995). The goal of the analy- 
sis is the decomposition of the observed line-of-sight ve- 
locities into inward and outward propagating waves in an 
annular region centered on the flux tubes. This allows us 
to detect the effect of a magnetic feature on the wave field 
through the difference between the outward and inward 
radially propagating Hankel components. In a spherical 
polar coordinate system (8, <f>) the wave components take 
the form 

¥ m (0,^,t) = e *(*»*+*«*) x 
x[A m {L,v)H£>{L0) + B m {L,v)H${L0)\, (1) 

where m is the polar azimuthal order, Hm and ffi are 
Hankel functions of the first and second kind respectively, 
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t is time, v is the temporal frequency, L = [1(1 + 1)] 1//2 
where I is the spherical harmonic degree of the mode, 
and A m and B m are the complex amplitudes of incoming 
and outgoing waves respectively. For small holds of view, 
such as the data we consider in this paper, the expan- 
sion is well approximated in terms of Hankel functions 
as opposed to Legendre functions. 

The Hankel analysis was centered on the locations of 
small magnetic elements as well larger set of con- 
trol locations offset by fixed distances from each mag- 
netic element. There were no active regions present in 
the MDI observations used. We hereafter refer to the 
control locations as "quiet-Sun" regions, although no at- 
tempt was made to select (or exclude) these locations 
based on magnetic properties of the MDI pixels. Thus 
the primary difference between the "magnetic feature" 
and "quiet Sun" locations is that the former coincide 
with the peaks of small magnetic flux regions while the 
latter are pseudo-random locations. The data consist of 
the set of 102 "high-resolution" MDI Doppler and mag- 
netogram datacubes, each with a 4 hr i nterva l, which 
has been previously used by IDuvall et al.l (|2006[ ) to mea- 
sure travel-time kernels for time-distance helioseismol- 
ogy. Each of the 102 regions are confined to the MDI 
"high resolution field" which spans 11 x 11 arc minutes 
and is centered about 160 arc seconds north of disk center 
(|Scherrer et al.llT995[ ). The data are obtained from all of 
the MDI full-resolution observations of at least 4 hours 
in duration and spanning the years 1996 and 1997. The 
cadence of both the Dopplergrams used in the helioseis- 
mic analysis, and the magnetograms used for identifying 
small magnetic elements, is 60 seconds. 

The criteria and procedure for the selectio n of the mag- 
netic elements are described in detail by IDuvall et al.l 
IpOOfih . The average magnetic feature has a peak mag- 
netic flux density of 76 G, and a full-width-at-half- 
maximum (FWHM) of 2.6 Mm, as determined from a 
fit by a 2D Gaussian funct ion to the averaged magnc- 
togram (|Duvall et al.l 120061 ). We used locations of the 
features from tables provided to us by Duvall. In to- 
tal, nearly 3400 locations of features were used. This 
is a somewhat larger sample than the (a pproximately) 
2500 features used by IDuvall etHI (|2006l) since we used 
features closer (i.e. as close as 40 pixels) to the edge 
of the datacubes than used in that work. There were 
four quiet-Sun locations offset from each magnetic fea- 
ture used, thus yielding a quiet-Sun control sample of 
13,600. Figure Q] shows the absolute value of the time- 
averaged magnetic flux density corresponding to one of 
the Doppler datacubes, showing the sample of features 
for that datacube. 

The Fourier-Hankcl decomposition was performed for 
all locations, and the coefficients A m (L, v) and B m (L, v) 
were determined for waves within an annular domain 
with colatitude relative to the center point of the analy- 
sis 6 ranging from 9 m i n = 10 pixels to Q m m = 40 pixels, 
where a pixel corresponds to 0.034 heliocentric degrees 
(or 0.413 Mm). Detail s of th e decomposition method 
are described by IBraunl (|1995| ). The method consists of 
discrete numerical transforms in the azimuthal, colati- 
tude and temporal domains. The azimuthal transform 
is computed for integer values of azimuthal order -10 
< m < 10 at the highest values of L. Not all of these 
orders m are computed (or useful) for smaller values of 
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Fig. 1. — The absolute value of the magnetic flux density, after 
averaging over 4 hours, for one of the 102 MDI regions. White 
circles indicate the position of the magnetic elements used for the 
Hankel analysis and their size corresponds to the outer annulus. 

wavenumbers fsee IBraunlll995l for details). The values 
of L and v for which the relative colatitude and tempo- 
ral transforms are computed compose a grid with spacing 
AL = 2Tr/(e max -e min ) = 352.9 and Av = l/T = 0.0694 
mHz, where T is the duration of the observations (4 hr). 

We focus in this paper on comparisons between obser- 
vations and models of phase-shifts between the outgoing 
and incoming wave components. Our primary motiva- 
tion for this is that the measurement and interpretation 
of observations of amplitude differences (e.g. due to ab- 
sorption) arc affected by factors such as the presence of 
background convective (or instrumental) noise and by de- 
tails of the solar excitation and damping mechanisms of 
the waves. We make no attempt to realistically include 
these in the numerical models described later. This chal- 
lenges our ability to make meaningful comparisons of ab- 
sorption coefficients, for example. Some of these effects 
further restrict the measurement of phase-shifts to those 
modes (typically with lower wavenumbers and tempo- 
ral frequencies) which have lifetimes significantly longer 
than the time the waves ta ke to propagate across the 
entire annulus (|Braunlll995f ). 

To determine feature-averaged phase shifts, we com- 
pute the summations, over the magnetic element and 
quiet Sun ensembles, of the product B m (L,v)A t n . l (L,v) 
where the asterisk denotes the complex conjugate. The 
effects of a "spurious phase-shift' caused by the leak- 
age of wave a mplitudes acr oss the wavenumber-frequency 
domain (see IBraunl I1995D are removed by considering 
the relative difference in the phase-shift between the 
magnetic-element and quiet-Sun ensembles. Thus, we 
compute a "corrected" ensemble- averaged product, 

= l(^n^) m Je i; ^(^A*J m J-arg((i 3m A;j qa )] ; 

where the brackets indicate the ensemble average over 
the samples of magnetic elements (me) and quiet-Sun 
(qs) locations and the explicit dependencies of A m and 
B m on L and v are omitted for clarity. 

An additional averaging across the width of the /-mode 
ridge, at each wavenumber L, is performed, such that the 
resulting phase shift is given by 

/ /-l/o (£)+*!/ \ 

S ohB (L) = arg / (B m (L, v)A* m (L, v))' me dv 

\Jv (L)-8v ) 

(3) 

where v$(U) is the /-mode frequency at wavenumber L. 
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The window 8v is determined by an inspection of the 
power spectra \A m (L, v)\ 2 such that it contains most of 
the /-mode ridge power. The value of 5v increases from 
0.3 mHz at L = 706 to 0.7 mHz at L = 1765. Figure 
[2] shows the observational phase shift, including its vari- 
ation with L for several values of m, and its variation 
with m dX L = 1412. In the next section, we describe the 
model we construct to reproduce these phase shifts. 
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Fig. 2. — Top panel: Variation of the observational phase shift 
with L for the azimuthal order m = (solid line), m = 1 (dashed 
line), and m = 2 (dashed-dottcd line). Bottom panel: Variation 
of the observational phase shift with m at L = 1412. The error 
estimates are obtained from an average over m of the absolute value 
of the phase-shift difference between +m and — m, divided by \/2. 
This estimation assumes that the error does not depend on m. 



3. NUMERICAL PROCEDURES 

We have solved numerically the nonlinear thrcc- 
dimensional (3D) MHD equat i ons using the code M ancha 
HKhomenko fc Colladosl [20061: iFelipe et alJl2010h . The 
code solves the equations for perturbations, which are 
obtained by subtracting the equations of initial magne- 
tohydrostatic equilibrium from the system of nonlinear 
MHD equations. Spatial derivatives are discretizcd using 
fourth-order centered differences and the solution is ad- 
vanced in time using a fourth-order Runge-Kutta scheme. 
It is stabilized by artificial diffusion terms and its paral- 
lel design is based on a domain decomposition scheme. 
A perfectly ma tched layer (PML) boundary condition 
(|Berengerlll994l) is applied in order to avoid wave reflec- 
tions. 

As a magnctostatic background, we ha ve used flux tube 
model s constructed using the method of [P neuman et al 
following the routines by IKh 
We consider 
tube 



1986), 



2008) 



IKhomenko et al 



routines by 
three flux tube models: 
radius of 170 km and 



the 
pho- 



small flux tube has a radius ol 1YU km and a 
tospheric magnetic field strength around 1600 G, with 
slight variations with height and radial distance inside 
the tube; the medium flux tube has the same mag- 
netic field strength, but a radius of 370 km; and the 
larger flux tube has also 1600 G photospheric field 
strength and a radius of 560 km. The external atmo- 
sp here outside the flux tube is the quiet Sun model 
S ()Christensen-Dalsgaard et al.lll996|) stabilized follow- 
ing the method described by iParchevskv fc Kosovichevl 



(|2007l) to avoid convective instability. The radial transi- 
tion between the magnetic and non-magnetic atmosphere 
is performed smoothly using a cosine which reduces the 
magnetic field from its maximum to zero over 100 km in 
order to avoid numerical problems due to the disconti- 
nuity in the magnetic field. For the small flux tube the 
magnetic field is strictly zero for radial distances higher 
than 250 km, in the case of the medium tube it vanishes 
at 450 km, while for the larger tube the magnetic field 
is zero beyond a radius of 630 km. Figure |3] shows the 
characteristic velocities, the ft = P ga s I Pmag parameter, 
where P ga s is the gas pressure and P ma g the magnetic 
pressure, and the pressure scale height at the axis of the 
tubes. 
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Fig. 3. — Properties of the flux tubes at the axis. Top panel: 
sound speed (solid line) and Alfven speed (dashed line) for the 
560 km radius tube; middle panel: /3 parameter for the tubes with 
radius of 170 km (dashed line), 370 km (solid line), and 560 km 
(dotted line); bottom panel: pressure scale height. 

We use a local Cartesian geometry defined by the hor- 
izontal coordinates x and y and the vertical coordinate 
z. The computational domain spans from z = —6 Mm 
to z = 0.6 Mm, where z = corresponds to the height 
where the optical depth is unity at a wavelength of 5000 
A in the quiet Sun atmosphere. The horizontal extent of 
the domain is a; € [-42.3, 32.3] Mm and y e [-32.5, 32.5] 
Mm, with the axis of the vertical flux tube located at 
x = 0, y = Mm. The spatial step is 50 km in the three 
spatial dimensions. In the vertical direction, PML layers 
with a thickness of 5 grid cells were set in the top and the 
bottom boundaries, so the physical domain spans from 
z = —5.75 Mm to z = 0.35 Mm. In the x direction a 
10 points PML was used, while in the y direction we set 
periodic boundary conditions with no PML. 

We aim to study the scattering of an /—mode by the 
flux tube. The vertical velocity of an /—mode wave 
packet which propagates from left to right in the x di- 
rection in a horizontally homoge neous atmosphere is de- 
scribed bv lCameron et all (|2008D as: 



v z (x,y,z 



t) = ReJ2 A ke kz e- ik{x ~ X0 



)-iu> k t 



(4) 



where Re means real part, are complex amplitudes, 
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Xq indicates the initial position of the wave packet, t is 
time, and ujk = \fgo~k is the eigenfrequency at wavenum- 
ber k. This wave packet is uniquely determined by the 
initial conditions: 

v = ReY2{i* + i)A k e kz+lk{x - X0 \ (5) 
k 

Pl = ReJ2i^ 1 Po9oe kz+ik{x - xo) , (6) 

k 

pi = ReY,iA k Uk 1 PoH p - 1 e kz+ik ( x - x °\ (7) 
k 

where pi and pi are the Eulcrian perturbations in pres- 
sure and density, respectively, po is the unperturbed den- 
sity, g = — goz is the acceleration due to gravity, z is a 
unit vector pointing upward, and H p is the density scale 
height given by the expression: 



1 



1 



Tq dz 



(8) 



where To is the background temperature, H p = c 2 s /(^go) 
is the pressure scale height, and 7 is the ratio of specific 
heats, and c s is the sound speed. Equations and 
(J7J arc obtained from the initial d isplacement vector for 
an /—mode (Equation (19) from iCameron et al.l 12008ft 
and its relations with the perturbations in density and 
press ure {e.g., Equations (4) and (5) from lCameron et al.1 
120081 respectively). 

Following Eqs. JS])-© we have introduced an /—mode 
at t = s located at xq = —37.3 Mm and spanning 
along the full domain in the y direction. As an initial 
distribution of /—mode amplitudes as a function of 
L = kR@ , where k is the horizontal wavenumber and Rq 
is the solar radius, we have imposed a Gaussian centered 
at L = 1000 and with a half width of 600. This is an 
unrealistic distribution, but we are interested in quan- 
tifying the absorption and phase shift at each L rather 
than reproducing the solar spectrum. We have chosen 
the distribution of /—mode amplitudes in order to get 
enough power in the wavenumbers of interest. Since we 
took real values for A k , all waves with different k are 
in phase at the starting position. We have imposed a 
small amplitude in order to be sure that the simulation 
remains in the linear regime. The initial location of the 
wave packet xq was selected to introduce the wave out- 
side of the outer circumference of the Hankel analysis. 
The duration of the simulation is T = 180 min, which 
corresponds to the total time that the wave packet needs 
to travel through all the domain in the x direction. The 
output is saved every minute, producing a set of 181 
three-dimensional cubes that provide the temporal evo- 
lution of the three components of the velocity, pressure, 
density, and the three components of the magnetic field 
for all the computational domain. 

Together with the 3D simulations of the atmosphere 
with the flux tube, we have also performed a two- 
dimensional simulation in a quiet Sun atmosphere, with- 
out the flux tube being present, but otherwise using ex- 
actly the same configuration as in the flux tube model 
computation. According to the set up of the 3D simu- 
lation, all the grid rows in the x direction for a corre- 
sponding y arc equivalent, except for the presence of the 



tube, which allows us to use the 2D quiet Sun simulation 
as a reference to obtain the scattered wave as the dif- 
ference between both simulations by subtracting the 2D 
simulation from all the xz planes in the 3D computation. 

4. F— MODE SCATTERING 

When the /—mode reaches the flux tube different tube 
waves are excited. These waves propagate upward and 
downward along the magnetic field lines, extracting en- 
ergy from the acoustic cavity and producing an effective 
absorption of the /—mode energy. The tube waves are 
visible in Figure^ which shows the z velocity (left panel) 

1/2 

and x velocity (right panel) scaled with factor p at 
t = 100 min in a region in and around the flux tube. As 
these waves propagate downward their amplitude drops 
due to the higher density at deeper layers, while their 
wavelength also decreases because of the lower Alfven 
velocity. 
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Fig. 4. — Vertical cuts of the ^-component of velocity (left panel) 
and the x-componcnt of the velocity (right panel) scaled with factor 
1 /2 

p„ at i = 100 min in the presence of a 370 km radius flux tube 
after substracting the quiet Sun simulation. The white-black colors 
mean positive-negative velocity directions; the range of the grey 
scale is the same in both panels. Vertical solid lines represent the 
boundaries of the flux tube. Vertical dashed lines correspond to 
the outer limits of the jacket modes region (see FigurcO, although 
they are only visible close to the tube. 

Figure [5] shows the vertical velocity at two different 
time steps for the simulation with the medium flux tube 
with R — 370 km. The left hand column corresponds 
to t = 80 min, just before the main part of the wave 
packet reaches the flux tube, while the right hand column 
illustrates the simulation at t = 100 min. The top panels 
represent the scattered waves in a xy plane at z = —0.5 
Mm, while the middle and bottom panels show the full 
wave field in a xy cut at the same height and a xz cut at 
y = Mm, respectively. Note that, in order to improve 
visualization, the color scale is different in the top two 
panels, which are 10 times more saturated. 

The scattered wave is obtained by subtracting the quiet 
Sun two dimensional simulation from every xz plane 
along the y direction of the flux tube simulation. In 
the 170 km radius flux tube the amplitude of the vertical 
velocity of this wave is 0.020 times the amplitude of the 
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incident /—mode, while for the medium tube (with 370 
km radius) the ratio between both amplitudes is 0.093. 
In the case of the big tube (with 560 km radius) the ra- 
tio of the amplitudes is 0.235. Since all tubes have the 
same magnetic field strength, their magnetic flux is pro- 
portional to R 2 , where R is the tube radius. We find 
that the scattering amplitude approximately scales with 
the magnetic flux, that is, with R 2 . 

The oscillations of the waves which travel along the 
tube are basicaly the generators of the scattered wave, 
which mainly corresponds to an /—mode, although a 
small amount of power is also scattered into the pi— mode 
(see Section [6.3|1 . Since the scatterer is axisymmctric, no 
scattering can be produced from an azimuthal order m to 
a different one. In Figure [5] it is seen that the scattered 
wave produced by the 370 km radius tube is composed 
by a mixture of m = and m = ±1 waves, dominated 
by the later ones. 

5. JACKET MODES 

In addition to the excitation of tube waves and the 
scattering into different wave modes presented in the 
previous sections, a scattering object can also generate 
jacket modes in its surroundings. Jacket modes were first 
discussed in the c ontex t of solar acoustic oscillations by 
IBogdan fc Calrvl (fl995h . They consist of a continuous 
spectrum of horizontally evanescent wave modes which 
propagate vertically in the non-magnetic region around 
the flux tube. The jacket modes are necessary because 
the vertically evanescent /- and p-modes alone cannot 
ensure the continuity of pressure and horizontal velocity 
of the oscillations across the flux tube boundary, due to 
the presence of the tube waves that propagate downward 
at large depth. 

The jacket modes obtained in the medium tube simu- 
lation can be seen in Figure|3]as the downward propagat- 
ing waves outside the boundaries of the tube with small 
vertical wavelength, especially in the horizontal velocity 
(right panel). 

The energy extracted by the flux tube from the 
/—mode goes to the tube waves and the jacket modes, 
which transport the energy upward and downward re- 
moving energy from the acoustic cavity. It is interesting 
to evaluate how the energy is distributed among these 
modes. In this analysis we will only consider the en- 
ergy which goes to deeper layers, represented by a nega- 
tive energy flux, since the proximity of the top boundary 
hinders obtaining a reliable positive flux into the atmo- 
sphere. The wave energy fl uxes were calculated following 
IBrav fc Loughheadl (|1974l) . The acoustic energy flux is 
given by the expression: 

F ac = piv, (9) 

and magnetic wave energy flux is given by: 

■ mag 

= Bix(yx B )M). (10) 

In these expressions v and Bi are the perturbed velocity 
and magnetic field, respectively, Bo is the background 
magnetic field and fig is the magnetic permeability. 

In the case of the waves inside the tube, the acoustic 
flux inside the tube can be obtained as the difference 
between the acoustic flux in the simulation with the flux 
tube being present and the quiet Sun simulation. With 
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Fig. 5. — Vertical velocity at t = 80 min (left panels) and t = 100 
min (right panels) for the simulation with a 370 km radius magnetic 
flux tube. Top: horizontal cut at z = —0.5 Mm of the scattered 
wave; middle: horizontal cut at z = —0.5 Mm of the full wave 
field; bottom: vertical cut at y = Mm for the full wave field. The 
white-black colors mean positive-negative velocity directions; the 
grey scale in the top panels is 10 times more saturated. In the top 
and middle panels the circle at x = Mm and y = Mm indicates 
the location of the tube, while the dashed line is the position of 
the cut shown in bottom panels. In the bottom panels, the vertical 
solid lines are the boundaries of the tube and the horizontal dashed 
line is the position of the cuts shown in top and middle panels. In 
the full wave field plots the tube modes are hardly visible because 
their amplitude is small compared to the /—mode. See Figure [4] 
for a dctalied plot of the tube wave velocities. 

regards to the magnetic wave flux, it is directly obtained 
from the application of Equation (|10[) to the region where 
the magnetic field is different from zero. Adding both 
fluxes we retrieve the total energy flux carried downward 
by the tube waves. At each height, we have summed the 
energy flux correponding to all the points inside the tube 
for all the time steps between the time that the /—mode 
wave packet reaches the tube and the time that the wave 
packet leaves it. The result is shown in Figure [51 where 
the variation of the wave energy flux of the tube waves 
with the height is plotted. Around z = km the total 
flux vanishes. Above that height the flux is positive, 
which means that the energy propagates upward, and 
below z — Mm it is negative, showing a minimum 
around z = —0.2 Mm in the case of the two smaller 
tubes and around z = —0.5 Mm for the larger tube. As 
these waves propagate downward they are damped by the 
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diffusivity, and their energy flux tends to zero at about 
z = —2 Mm in the case of the 370 km and 560 km radius 
tube and at z = — f Mm for the 170 km radius tube. The 
downward wave energy flux for the two bigger tubes is 
similar, being 1.6 times higher than the corresponding to 
the small tube. 




-1.5 -1.0 
HEIGHT [Mm] 

Fig. 6. — Total vertical wave energy flux inside the magnetic flux 
tube with 560 km radius (dotted line), the tube with 370 km radius 
(solid line), and the tube with 170 km radius (dashed line). All 
cases are normalized to the absolute value of the flux of the 560 
km radius tube at z = —0.5 Mm. 

The jacket modes appear in the non-magnetic region, 
surrounding the magnetic flux tube. Thus, only the 
acoustic wave flux contributes to its wave energy. How- 
ever, some care must be taken to obtain a correct evalu- 
ation of its energy flux. Outside the flux tube, the scat- 
tered wave field (obtained as the difference between the 
simulation with the flux tube and the quiet Sun simula- 
tion) is composed by the scattered wave (fundamentally 
an /—mode) and the jacket modes. The difference be- 
tween the acoustic flux in both simulations will corre- 
spond to the sum of the scattered wave and the jacket 
mode fluxes. To isolate the contribution of the jacket 
modes, we have filtered out the low vertical wavenum- 
bers of the scattered velocity in the surroundings of the 
flux tube, and we have calculated its acoustic flux using 
that velocity in Equation |U The filter selects the waves 
with vertical wavenumber higher than 3.14 Mm -1 and, 
thus, only waves with vertical wavelength below 2 Mm 
are considered for the jacket modes. At the photosphere 
the jacket modes have a vertical wavelength around 1 
Mm. As they propagate downward, their wavelength is 
reduced in order to match the decrease of the wavelength 
of the tube waves due to the reduction of the Alfven 
speed, showing a 0.45 Mm vertical wavelength at z = —3 
Mm. 

We have summed the energy flux of the jacket modes in 
an annular region surrounding the flux tube at the height 
where the downward tube wave energy flux is maximum, 
that is, at z = —0.2 Mm for the two smaller tubes and 
at z = —0.5 Mm for the larger tube. The annulus is 
delimited by the radius of the flux tube in the inner part, 
and the size of the outer radius varies. Figure [7] shows 
the result. In the case of the 370 km radius flux tube 
(solid line) the energy flux of the jacket modes increases 
with the outer annulus radius until R = 1 Mm. For larger 
radius the energy is constant, which means that at larger 
distances from the tube there is no energy flux associated 
with the jacket modes. Thus, the size of the jacket mode 
region is about 0.6 Mm around the tube. The energy flux 



has been normalized to the absolute value of the energy 
of the waves inside the tube at the same height. We find 
that the energy carried downward by the jacket modes is 
approximately 40% of the tube waves energy for the 370 
km tube. With regards to the 170 km and 560 km flux 
tubes, a similar size of the jacket mode region is obtained. 
However, in these cases the energy transported by the 
jacket modes is around 16% of the energy extracted by 
the waves inside the tube. 
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Fig. 7. — Total vertical acoustic wave energy flux of the jacket 
modes inside an annular region surrounding the flux tube at z = 
—0.2 Mm. The inner radius of the annulus is given by the radius 
of the flux tube and the outer radius by the abscissa coordinate. 
Dotted line: 560 km radius; solid line: 370 km radius; dashed line: 
170 km radius. The two later cases are normalized to the energy 
flux of the tube waves at z = —0.2 Mm, while the 560 km case is 
normalized to the energy flux of the tube waves at z = —0.5 Mm. 



6. HANKEL ANALYSIS 

The numerical simulations were also analyzed using 
the Fourier-Hankel spectral decomposition method, as 
explained in Section [2] The coefficients A m {L,v) and 
B m {L,v) are evaluated from the vertical velocity at 
z = 0.2 Mm. A different sampling in L and v was ob- 
tained because of the differences in the annular domain 
and duration of the temporal series used. In this case, the 
analysis is restricted to the annular region delimited by 
the polar angles m in — 0.00589 rad and 9 max = 0.04310 
rad, which correspond to a radial distance of R m in = 4.1 
Mm and R ma x = 30 Mm, respectively. Thus, we obtain 
a grid with spacings AL = 2tt /(9 max — 6 m i n ) = 168.8, 
while the duration of the simulations T = 180 min pro- 
vides a Av = l/T = 0.0936 mHz. The outer radius of 
the annulus is given by the horizontal size of the compu- 
tational domain. It was chosen as a compromise between 
a big enough domain to obtain good sampling in L, but 
not too big to avoid very computationally expensive sim- 
ulations. 

To determine the absorption coefficient, the power of 
the ingoing and outgoing Hankel components has been 
summed across the ridge of the /—mode in order to re- 
trieve a better signal-to-noise level. At each given L, the 
power for the ingoing wave is determined as 

/■i/o(£)+<5f 

\A m {L)\ 2 = / \A m (L,v)\ 2 dv. (11) 

J vq (L) — 5v 

The frequency window Su was set to approximately 0.3 
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mHz, which corresponds to 3-4 frequency bins. The same 
average is applied to the outgoing B m (L, v) components. 

The absorption coefficient a m (L) along the ridge of the 
/—mode is then obtained as 



V J 



1 



\B m {L)\ 2 /\A m {Lf 



(12) 



and the phase shift is given by Equation (|3|), without 
applying the ensemble sum over magnetic elements (me) . 

The definition of the absorption coefficient may not 
correspond to a real dissipation of energy if there is sig- 
nificant mode mixing present. 

We have calculated the absorption coefficient of the 
simulation with the flux tube (a^ r (L)) as well as the 
quiet Sun reference simulation (a^ s (L)). Although in 
the latter one the absorption should vanish, since there is 
not any scattering element, due to numerical reasons (nu- 
merical diffusivity and interaction with the top bound- 
ary) some absorption appears at high wavenumbers. The 
absorption coefficient measured directly from the quiet 
Sun simulation is given by 



„os fn _ i _ \ B m ( I/ )l 2 fT /n 
m() \A% s (LW a{L) > 



(13) 



where A^f(L) and B^ S (L) represent the ingoing and 
outgoing power is the quiet Sun simulation, respectively, 
in the ideal case where there is no numerical damping. 
The numerical damping which produces a different power 
in the outgoing waves is included in cr(L). In the same 
way, from the flux tube simulation we measure the ab- 
sorption coefficient as 



FT 



(L) 



a(L). 



(14) 



We assume that the numerical damping a(L) is the 
same in both simulations, since they use exactly the same 
configuration. In the quiet Sun the ingoing power is equal 
to the outgoing power (\A^ S (L)\ 2 = \B$°(L)\ 2 ). Thus, 
from Equation (fT3|) we retrieve cr(L) = 1 — a^ s (L). Tak- 
ing into account that by definition the real absorption 
coefficient produced by the tube is 



m = i 



(15) 



including Equation (fT4")l and the expresion for a(L) ob- 
tained from Equation (|13[) in the previous equation, we 
obtain 



..FT 



a m (L) 



(L) 



1 



Q 



QS 



(L) 



(16) 



In the following we will discuss the absorption coeffi- 
cient obtained after applying this correction. 

In Section [4] we discussed qualitatively the properties 
of the scattered wave. We have performed the Hankel 
analysis of that wave by decomposing in Hankel functions 
the difference in the photosphcric vertical velocity at z = 
0.2 Mm between the wave field of the simulation with the 
flux tube being present and the quiet Sun simulation. In 
the simulation with the 170 km radius magnetic flux tube 
the analysis reveals that the power of the components 



with azimuthal order m = ±1 is 1.70 times higher than 
the power of the axisymmetric components with m = 
0. The power in higher azimuthal orders is negligible. 
For the medium tube (with radius of 370 km) the power 
of the m — ±1 components is only 1.19 times higher 
than the corrcponding to m — 0. This simulation shows 
a small amount of power in the azimuthal orders m = 
±2. On the other hand, in the 560 radius flux tube the 
dominant azimuthal order of the scattered wave is the 
axisymmetric m = 0, whose power is 1.1 times higher 
than the power of m = ±1. It shows a significant amount 
of power inm = ±2, which is 7.43 times smaller than the 
power in ttl = 0. 

6.1. f—mode absorption 

The top panel of Figure [8] shows the absorption coef- 
ficient as a function of L for the two lowest azimuthal 
orders determined for each of the magnetic flux tubes. 
On the one hand, in the case of the flux tube with a 
radius of 170 km (asterisks) the highest absorption is 
retrieved for the azimuthal order m = ±1. It reaches a 
value above 0.1 at L = 2500. The azimuthal order m = 
also shows a significant amount, presenting an absorption 
higher than half of the corresponding to m = ±1. For 
the rest of the azimuthal orders the absorption is very 
low, although high wavenumbers show some absorption 
at m = ±2. The variation of the absorption with the 
azimuthal order is clearly shown in the bottom panel of 
Figure [8] for L = 2364. It exhibits a perfect symmetry 
around m = 0, with clear peaks at m — ±1 and dropping 
to zero at higher m. 

On the other hand, for the medium flux tube with a 
radius of 370 km (diamonds) the measured absorption 
coefficient is generally higher. The highest absorption is 
also obtained for \m\ = 1, but in this case its value is 
around 0.3 at L = 2500. Below L = 1500, the azimuthal 
orders m — and m = 2 (not shown in the figure) pro- 
duce a similar absorption, but above that L they split up. 
The absorption of the later one keeps increasing with L 
and its a is around 0.1 at L = 2500, while the increase of 
the absorption in m = seems to be slower and it shows 
a maximum absorption around 0.05 at L = 2500. The 
variation of a with m (bottom panel of Figure [5]) also 
shows symmetry around m = 0. 

With regards to the 560 km radius flux tube, the high- 
est absorption is also retrieved for m = 1 , which presents 
an a = 0.47 at around L = 2500. At high L values 
the absorption coefficient shows significant absorption at 
m = 2 and even at m — 3, opposite to the m = case, 
which presents a very low a coefficient. 

Several conclusions can be extracted from the compar- 
ison of the absorption coefficient measured for different 
flux tubes. Firstly, the variation with L shows a similar 
tendency for all flux tubes. It increases almost linearly 
with L for all the azimuthal orders. Secondly, the ab- 
sorption is subject to the magnetic flux of the scattering 
element. A higher magnetic flux produces higher ab- 
sorption, although its distribution in wavenumber and 
azimuthal order depends on the radius of the scatterer. 
Despite the fact that the big flux tube has almost 11 
times higher magnetic flux than the smaller one, its ab- 
sorption coefficient in m — 1 is far from being infered 
as 11 times higher than the m = 1 absorption of the 
small tube, since the ratio between the absorption of 
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Fig. 8. — Top panel: Variation of the absorption coefficient with 
L for the azimuthal orders m = (dashed line) and m = 1 (solid 
line). Bottom panel: Variation of the absorption coefficient with m 
at L = 2364. In both panels asterisks correspond to the simulation 
with a flux tube with a 170 km radius, diamonds to the flux tube 
with a 370 km radius, and crosses to the flux tube with 560 radius. 

both tubes varies with L. At L = 1013 the to = 1 
absorption of the big flux tube is 16.5 times higher than 
the corresponding to the small tube, while at L = 2532 
the ratio is just 4.3. Finally, the relation between the 
absorption at different azimuthal orders is different. As 
can be seen in the bottom panel of Figure HI the higher 
the radius of the tube the higher absorption coefficient 
in all azimuthal orders except m = 0. Surprisingly, a 
in m = decreases with increasing radius, despite the 
higher power of the scattered wave in to — for the 
larger radius tubes discussed in the previous section. 

6.2. Phase shifts 

We arc interested in the difference between the ingo- 
ing and outgoing pha ses produced by th e scattering el- 
ement. According to iBraun et al.l (fl992). for a reliable 
determination of the phase shift between the incoming 
and outgoing waves it is necessary that the observations 
last long enough so that the wave packet can travel a dis- 
tance comparable to the annulus diameter. The temporal 
duration of our simulations (T = 180 min) was chosen 
in order to satisfy this condition. The phase shift was 
evaluated following Equation ([3|). The values discussed 
in this section correspond to the difference between the 
simulation with the flux tube and the quiet Sun reference 
simulation. 

The top panel of Figure [9] shows the variation of the 
phase shift with L for the azimuthal orders m = 1,2. 
Starting with the smaller tube (asterisks), we find that 
the phase shift increases from 8 m= i = 0° at L = 800 
to S m= i sa 10° at L = 2500. The variation of the az- 
imuthal order to = with L (not shown in the plot) is 
very similar to to = 1 , showing an almost linear increase 
with an slightly lower phase shift. Below L = 800 the 



phase shifts of both azimuthal orders are around 0. The 
azimuthal order m = 2 shows a much smaller phase shift, 
and it almost vanishes for all L values. As in the case 
of the absortion coefficient, significant phase shifts are 
only obtained for m = and to = ±1 azimuthal orders 
(Figure [H bottom panel). 

The behavior of the phase shift produced by the 
medium flux tube (diamonds) is similar to the smaller 
one, but showing a much higher value. In this case, the 
phase shift is also around 0° below L = 500 and it in- 
creases with L until reach 8 ~ 35° for both m = and 
m = 1 azimuthal orders. As in the case of the small flux 
tube, the phase shift of the to = 1 azimuthal order is 
a bit higher, and this difference increases with L. The 
phase shift corresponding to m = 2 is much lower, since 
it only reaches 8 m= 2 ~ 7°. 

Finally, the azimuthal order to = 1 of the larger tube 
(crosses) shows an increasing phase shift which reaches 
almost 8 ps 70° at L = 2500. In this case, the phase shift 
in m = is slightly higher than for m = 1, as shown in 
the bottom panel of Figure [9l The phase shift produced 
in to = 2 is more significant than in the other tubes, 
since it is around 5 ss 28° at L = 2500. 
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Fig. 9. — Top panel: Variation of the phase shift with L for 
the azimuthal orders m = 1 (solid line) and m = 2 (dashed line). 
Bottom panel: Variation of the phase shift with m at L = 2364. 
In both panels asterisks correspond to the simulation with a flux 
tube with a 470 km radius, diamonds to the flux tube with a 370 
km radius, and crosses to the flux tube with a 560 radius. The 
dashed line in the bottom panel corresponds to the 560 km radius 
flux tube at L = 4765. 

A few points deserve special attention. Although the 
absorption coefficient of to = decreases with the radius 
of the tube (Figurc[5]), it is interesting to note that for the 
phase shift the azimuthal orders has different relevance. 
In the bottom panels of Figure O it is clearly seen that 
the phase shift produced in to = has a significant role. 
In fact, the to = phase shift is very close to the one 
obtained for m = 1, and even higher in the case of the 
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larger tube. In all tubes the phase shift shows a similar 
tendency, being significant for m = 1 and m = and 
small for m — 2, despite the clear differences which are 
present in the absorption coefficient. However, the larger 
radius (and higher magnetic flux) of the bigger tubes does 
affect the value obtained for the phase shift. At L = 844 
the phase shift obtained for m = 1 in the big flux tube is 
13 times higher than the corresponding to the small flux 
tube. This difference decreases with L, and at L = 2532 
the ratio between the phase shift of both tubes is around 
7. Note that the magnetic flux of the big tube is around 
11 times higher than that of the small tube. A similar 
behavior is obtained for the azimuthal order m = 0, and 
also including the medium tube in the comparison. The 
azimuthal order m = 2 shows a different pattern, since 
the ratio of the phase shift between a bigger and a smaller 
radius tubes increases with L. 

6.3. Mode mixing 

The absorption coefficient gives us a measure of the 
power lost by a certain mode, which corresponds to a 
frequency, wavenumber and azimuthal order. However, 
it does not necessarily mean that part of its energy has 
suffered a real absorption. The scattered wave may be a 
different wave mode, with different degree L and radial 
order n or azimuthal order m. For structures that are 
stationary in comparison to the typical wave period it 
is assumed that the outgoing wave must has the same 
frequency of the incident wave. In our simulations this 
condition is strictly satisfied, since we are using a magnc- 
tohydrostatic model. The change of the incident mode n 
to a different scattered mode n' at a fixed frequency pro- 
duced by a magnetic element is commonly called mode 
mixing. In addition, since these flux tube models are 
axisymmetric, no scattering can be produced from an 
azimuthal order m to a different order m' 

We have tried to measure the mode mixing produced 
by the three magnetic flux tube models. Since we are 
introducing as an initial condition the propagation of an 
/—mode, we can only estimate the scattering produced 
from this incident /—mode to higher order modes. How- 
ever, from these simulations we have only retrieved a 
significant amount of power in the p\ ridge, so we have 
evaluated the scattering from the /—mode to the p\. We 
have defined the following quantity 



a/_piO) 



|5 p ™i 2 



\b9 s Mi* 



(17) 



where |B^/>)| 2 is the power in the p\ ridge for the out- 
going component of the quiet Sun simulation, \B^{u)\ 2 
is the power in the p\ ridge for the outgoing component 
of the simulation with the flux tube, and \B^ S (v)\ 2 is the 
power in the /—mode ridge for the outgoing component 
of the quiet Sun simulation. All these expressions corre- 
spond to the sum of the power for the azimuthal orders 
m = —1,0, 1 and the degrees L over which the ridge of 
corresponding mode spans. The coefficient oif— p x repre- 
sents the ratio between the power in the pi-mode gen- 
erated by the flux tube and the power in the /—mode 
in the case without the tube being present at the same 
frequency. A negative value means that there has been 



emission to the p± — mode. In the numerator we have in- 
troduced the difference in power in the pi-mode between 
the quiet Sun and the flux tube simulations in order to 
correct from the small amount of energy which appears 
in this mode from the initial condition. We have intro- 
duced as initial condition the analytical expression of a 
propagating /—mode. However, due to the limited size 
of our computational domain and the presence of PML 
layers the solution of an /—mode in our domain is slightly 
different, and a very small amount of power goes to the 
Px ridge. 



0.000 



-0.005 



-0.010 



-0.015 



-0.020 




2.5 



3.0 



3.5 4.0 
Freq [mHz] 



4.5 



5.0 



Fig. 10. — Variation of the coefficient ctj_ p i (Equation 1171 ) with 
frequency. Crosses: 560 km radius flux tube; diamonds: 370 km 
radius flux tube; asterisks: 170 km radius flux tube. 



The ctf-pi coefficient obtained for all simulations is 
plotted in Figure [lOj No measurable amount of power 
was found in the p\ ridge below 2.6 mHz and above 4.6 
mHz. On the other and, the low power obtained between 
these frequencies makes it hard to retrieve a reliable mea- 
sure of the mode mixing. For the large flux tube we find 
that the power scattered to the p\ mode is below 2 x 10 -2 
times the power of the incident /—mode. The coefficient 
otf- p i seems to decrease with the frequency, although 
the result is very noisy. In the case of the medium tube, 
atf-pi is lower than 3 x 10~ 3 . The smaller tube show a 
much lower mode mixing, and its a/_ p i decreases from 
at 3.4 mHz to -5 x 10" 4 at 4.4 mHz. 

6.4. Comparison with observations 

Figure [2] represents the phase shift obtained from the 
observations. As shown by the bottom panel, it has a 
broad distribution in azimuthal order, exhibiting a sig- 
nificant phase shift between m = —6 and m — 5, with 
a peak at m — 0. The variation of the phase shift with 
m greatly differs from the one measured for the simu- 
lation, where in the two smaller radius cases we found 
that the phase shift is concentrated in the azimuthal or- 
ders m = —1,0, 1 with symmetric distribution, peaking 
at m = — 1 and m = 1, while the larger tube shows some 
significant phase shift in m = 2 (Figure ^ . 

We consider a scenario which could explain this mis- 
match. We have tried to mimic the error introduced in 
the observational Hankel decomposition by a displace- 
ment in the position of the annulus with respect to the 
magnetic element. This spatial translation can strongly 
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Fig. 11. — Top panel: Variation of the phase shift with L for the azimuthal orders m = (solid line) and m = 1 (dashed line). The bars 
above each L with data indicate the errors, the bottom one corresponds to the simulation and the top one is the error from the observations. 
Bottom panels: Variation of the phase shift with m at different L: from left to right L = 1056, L = 1412, and L = 1765. In all panels 
asterisks correspond to the simulation with a flux tube with a 560 km radius after averaging and diamonds to the observations. 



affect the incoming and outgoing power retrieved from 
the Hankel decomposition. As an example, consider a 
point source of waves which generates a concentric pat- 
tern of wavefronts that propagates radially outwards, like 
the one created by a pebble dropped into a pond. When 
the origin of the coordinate axis is located at the position 
where the pebble hits the water, the Fourier-Hankel anal- 
ysis of the wave will reveal a large Bq (outgoing axisym- 
metric wave) and zero A (no incoming power). However, 
at different positions the power will appear in different 
azimuthal orders. As stated above, the observed phase 
shift is retrieved from 3400 locations of small magnetic 
elements. If the resolution limitations and the movement 
of the magnetic element during the 4 hours of observa- 
tion produce a slight shift in the determination of the 
center of the magnetic element, the average of the phase 
shift obtained from the 3400 elements, all of them with 
different displacement, will produce a broadening of the 
distribution of the phase shift with the azimuthal order. 
We have imitated this limitation by averaging the phase 
shift retrieved from 3400 realizations of the Hankel anal- 
ysis in the simulation with the 560 km radius tube, with 
the center of the annulus shifted randomly around the 
axis of the tube in a Gaussian distribution with FWHM 
equal to that obtained from an estimation of the Point 
Spread Function for MDI high-resolution da ta, which 
corresponds to 1.14 Mm (jTarbell et al.lll997| ). In this 
analysis we have used the same annular region used in 
the observations in order to obtain the same sampling in 
L. 

Figure QT] shows the phase shift obtained from the ob- 
servations (diamonds) and the corresponding to the sim- 
ulations after applying the average (asterisks). The top 
panel shows that both the simulated and observational 
phase shifts increase with L for the two azimuthal or- 



ders plotted. At L = 1059 the observational phase shift 
is higher than the simulated one, but the rest of the L 
show good agreement. The error in the simulations is ob- 
tained using the same estimation as in the observations, 
although it is very small due to the symmetry retrieved 
between the positive and negative azimuthal orders. 

In the bottom panels we can see how the simulations 
resemble the observational broad distribution with m. 
As higher L is considered, significant phase shifts are 
obtained at higher azimuthal orders. At L = 1412, in the 
simulations the phase shift drop to zero at m = ±5, while 
the observations present a slightly broader distribution, 
with higher phase shifts at higher m. The dashed line 
in the bottom panel of Figure |9] shows the phase shift 
obtained for the 560 km radius tube at L = 1765, the 
same L value of the right bottom panel of Figure QTJ 
Note the broadening produced by the shifted average of 
the Hankel analysis. Interestingly, it also produces the 
effect of generating a clear phase shift peak at m = 0, 
which differs from the original data, where the azimuthal 
orders m = and m = 1 shows approximately the same 
shift. 

7. DISCUSSION AND CONCLUSIONS 

We have presented the analysis of the scattering pro- 
duced by magnetic flux tube models using 3D nu- 
merical simulation s. Previous attemps to model this 
phenomenon (e.g., iGordoyskvv fc Jain! 12007 ; Uain et al.l 



[2551 lHanasoge et al.ll2008t iHindman fc Jaidl20r2f )" have 
faced the problem by means of analytical treatments. 
These types of studies are the first steps toward the 
comprenhension of this issue and they provide a valu- 
able heritage to understand the wave interaction with 
magnetic media and confront it with the forthcoming 
observations and modeling. However, the simplifications 
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needed to carry out their development restrict their re- 
sults to some idealized cases. From this scope, the use of 
numerical simulations emerges naturally as the next step 
to address these questions in more general situations. 

In these simulations we have propagated 
an /—mode through a mod el S atmosphere 
(|Christensen-Dalsgaard et all 119961) stabilized against 
convective instabilities embedded with a flux tube 
model. In order to compare how some properties of 
the tubes affect the scattering three realizations were 
performed, using different flux tube models. All tubes 
have the same peak magnetic field strength, but they 
differ in the radius and, thus, in the magnetic flux. 

Our simulations show that the interaction of an 
/—mode with a flux tube excites tube waves. These 
waves propagate along the magnetic field lines and pro- 
duce a real absorption of the incident energy, since it 
is extracted from the acoustic cavity. The oscillations 
of the tube produced by these waves generate a scat- 
tered wave. It is composed of a mixture of axisymmetric 
(m = 0) and dipolar (m = ±1) modes, whose distribu- 
tion in frequency and azimuthal order depends on the 
radius of the flux tube. For thin flux tubes, the m = ±1 
dipolar oscillation dominates the tube wave, while ax- 
isymmetric oscillations (to = 0) become important for 
larger tubes . This result agrees w ith those previously 
obtained by iDaiffallah et all (|2011l) . 

We have quantified the absorption coefficient and 
phase shift produced by the three magnetic flux tube 
models. Based on the results discussed in the previous 
section, we draw the following conclusions. Firstly, the 
absorption increases with wavenumber (frequency) for all 
azimuthal orders and tube models. Secondly, the amount 
of absorption in general increases with the magnetic flux 
of the tube, although this increase depends significantly 
on the wavenumber and azimuthal order. Thirdly, the 
distribution of the absorption in azimuthal order depends 
on the radius of the tube. In all models the peak absorp- 
tion is obtained for m = 1. However, in the tube with 
170 km radius it is followed by m = 0, with a weak ab- 
sorption in to = 2, while in the 370 km radius tube the 
absorption in m = 2 is stronger than the corresponding 
to Tii = 0, and in the case of the 560 km radius tube the 
absorption in m = is especially low. The absorption of 
the axisymmetric m = order decreases with the radius. 

It is noticeable that the different behavior that the 
phase shift shows regarding the second and third points 
of the previous paragraph. From the simulations with 
the tubes we find a similar distribution in the phase shift 
produced in different azimuthal orders, although it seems 
to approximately scale with the magnetic flux, with some 
dependence on the L value and the azimuthal order. In 
this way, m = and \m\ = 1 show a very similar phase 
shift, the later slightly higher in the two smaller tubes 
and the opposite in the larger tube, while the phase shift 
produced in m = 2 is very small, except for the 560 km 
radius tube. 

In this work we are not only interested in modeling the 
scattering process, but also in applying this knowledge 
to interpret observations. A deeper understanding of the 
wave interaction with small magnetic scatterers can yield 
a basis to infer the properties of the scattering elements, 
even at scales smaller than the observational resolution. 
We have compared the numerical results for the phase 



shift with observations of an ensemble averaging of thou- 
sands of small magnetic elements. In order to perform 
an equivalent comparison, the phase shift obtained from 
3400 realizations of the Hankel analysis of the simula- 
tions with the 560 km radius tube with a small shift in 
the position of the annulus was also averaged. 

The phase shift produced by our larger tube model 
after averaging shows a good qualitative agreement with 
the observed phase shift. Since the phase shift scales with 
the magnetic flux of the scattering element, we may con- 
sider that the phase shift of the observed elements could 
be produced by fluxtubes with magnetic flux around the 
corresponding to tubes with 560 km radius and 1600 G. 
The current work suggests one possible solution for the 
properties of the tube model, although other combina- 
tions of radius and magnetic field strength might also 
work. This kind of measurement seems to be a promis- 
ing method to infer the characteristics of small magnetic 
networks elements. However, some caution must be con- 
sidered in their interpretation. Some of the observed 
magnetic elements used in this study show strong asym- 
metries (see Figure [1} . When the scattering element is 
not axisymmetric the scattering is not restricted to oc- 
cur from an azimuthal order to to the same order m, 
but the scatterred wave can correspond to a different or- 
der m! . These nonsymmetric elements could contribute 
to the broadening of the distribution of phase shift with 
m. On the other hand, in the observed magnetograms 
it is hard to find enough completely isolated magnetic 
elements. From the 3400 elements used in the analysis, 
many of them present other small magnetic features in- 
side the 16.5 Mm annulus in which the Hankel decompo- 
sition was performed. The observational analysis might 
be contaminated by the scattering produced by these 
other elements. 

In the comparison between the observations and the 
simulations, we have assumed that the observational 
shift in the center of the Hankel analysis is restricted 
to a Gaussian with the FWHM of the PSF from high- 
resolution MDI. If the proper motions of the flux tube 
have larger extension, the results could be affected. How- 
ever, the tubes are presumably moving at the time scale 
of the granulation, which is not so different from the 
time scale of the wave period. This could mean that the 
approximation of a stationary tube is not so good, and 
limits the capacity of this work to address this issue. It 
would be interesting to extend the analysis to moving 
flux tubes in the future. 

These results provide a warning to be cautious when 
interpreting ensemble averages of observational data. In 
observations like those presented in this paper, where 
the contribution of an individual magnetic feature is too 
low to get a reliable measurement, the average of sev- 
eral cases is a compulsory procedure to obtain a strong 
enough signal. However, the individual and unique char- 
acteristics of each element, together with the limitations 
to perform the analysis using exactly the same configura- 
tion, can lead to a result which may point to conclusions 
that do not reflect the real observed structure. In the 
particular case studied in this work, from the observa- 
tional broad distribution of the phase shift with to, one 
could assume that a big magnetic feature is necessary to 
produce that dependence with the azimuthal order. The 
analysis of the averaged simulated flux tube reveals that 
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a similar measure of the phase shift dependence can be 
retrieved from a very different magnetic element, mak- 
ing difficult infering an irrefutable conclusion about the 
nature of the observed elements. 
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